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We  perform  a  quantum  simulation  of  an  Ising  model  with  long-range  interactions  and  both  trans¬ 
verse  and  longitudinal  fields  -  the  most  general  external  field  configuration  -  in  a  system  of  6  to 
10  trapped  atomic  ions.  Quantum  fluctuations  at  zero  temperature  drive  the  ground  state  spin 
ordering  through  several  classical  first-order  phase  transitions  as  the  strength  of  the  longitudinal 
field  is  increased.  The  Hamiltonian  under  study  generates  a  fractal  Devil’s  staircase  structure  and 
maps  onto  a  large  number  of  many-body  and  energy-optimization  problems,  showing  how  quantum 
simulation  can  potentially  find  solutions  that  are  classically  intractable. 


Many-body  quantum  systems,  such  as  high-Xc  super¬ 
conductors  [1]  or  complex  configurations  of  interacting 
spins  [2-5],  are  often  difficult  to  describe  analytically  or 
calculate  numerically  on  account  of  the  exponential  scal¬ 
ing  of  the  Hilbert  space  with  the  system  size.  For  in¬ 
stance,  predicting  the  behavior  of  some  simple  quantum 
spin  models  is  currently  limited  to  just  a  few  dozen  spins 
[6] ,  while  solving  a  fully-connected  frustrated  Ising  model 
is  known  to  be  an  NP-complete  problem  [7].  Instead, 
quantum  simulators,  in  which  a  well-controlled  quantum 
system  is  used  to  simulate  a  system  of  interest  [8,  9], 
may  be  employed  to  calculate  ground  state  or  dynamical 
properties  of  a  Hamiltonian  that  would  otherwise  prove 
classically  intractable  or  non-integrable. 

Quantum  simulators  require  excellent  coherence  prop¬ 
erties  and  high-fidelity  readout  and  control.  To  date, 
simple  quantum  simulations  have  been  performed  in  a 
variety  of  systems,  using  neutral  atoms  in  optical  lat¬ 
tices  [10-12],  nuclear  magnetic  resonance  [13],  photons 
[14,  15],  and  trapped  atomic  ions  [16-25].  Early  trapped 
ion  simulations  of  many-body  physics  demonstrated  tun¬ 
able  spin-spin  interactions  [16,  17],  engineered  entangled 
spin  states  through  dissipation  [20] ,  and  observed  the  on¬ 
set  of  a  quantum  phase  transition  as  the  number  of  spins 
was  increased  from  2  to  9  [21].  Initial  studies  of  frustra¬ 
tion  within  the  context  of  the  antiferromagnetic  (AFM) 
Ising  model  [18]  were  later  enhanced  by  tuning  the  range 
of  interaction  and  the  degree  of  frustration  in  a  system 
of  up  to  16  spins  [24]. 

In  this  Letter,  we  report  the  first  quantum  simulation 
of  an  AFM  Ising  model  with  long-range  interactions  and 
both  transverse  and  longitudinal  magnetic  fields.  The 
addition  of  the  longitudinal  field  allows  for  the  applica¬ 
tion  of  an  arbitrary  magnetic  field  relative  to  the  Ising 
couplings  and  the  creation  and  observation  of  new  spin 
phases  at  zero  temperaure.  As  the  strength  of  the  longi¬ 
tudinal  field  is  increased  from  Bx  —  0,  the  system  passes 
through  N / 2  first-order  phase  transitions  driven  by  quan¬ 
tum  fluctuations,  where  N  is  the  number  of  spins.  Deter¬ 


mination  of  the  ground  state  spin  ordering  for  differing 
longitudinal  fields  reveals  a  Wigner-crystal  spin  struc¬ 
ture  [26],  maps  on  to  a  number  of  energy  minimization 
problems  [27,  28],  and  shows  hints  of  a  complete  Devil’s 
staircase  [29]  which  would  emerge  for  convex  long-range 
interactions  as  N  — >  oo. 

The  system  is  described  by  the  Hamiltonian 

u  E  +  4°  +  By{t)  E  4°  (!) 

i<j  i  i 

where  J-ij  gives  the  strength  of  the  Ising  coupling  be¬ 
tween  spins  i  and  j ,  Bx  is  the  magnitude  of  the  longitu¬ 
dinal  magnetic  field,  By(t )  is  a  time-dependent  transverse 
field,  and  cry  is  the  Pauli  spin  operator  for  the  ?’th  par¬ 
ticle  along  the  a  direction.  We  tune  J,;j  >  0  to  generate 
AFM  interactions  between  the  spins  that  fall  off  with  dis¬ 
tance.  At  t  —  0  the  spins  are  initialized  to  point  along  the 
direction  of  the  total  magnetic  field  B  =  Bxx  +  By(0)y, 
with  By{ 0)  J,  which  is  the  instantaneous  ground  state 
of  the  Hamiltonian  in  Eqn.  1  to  good  approximation. 
The  system  will  remain  in  the  ground  state  if  the  trans¬ 
verse  field  By{t)  is  ramped  down  adiabatically.  The  re¬ 
sulting  spin  order  at  the  end  of  the  ramp  when  By  =  0 
would  then  reveal  the  ground  state  of  an  Ising  Hamil¬ 
tonian  with  long-range  interactions  and  a  longitudinal 
field. 

Fig.  1(a)  shows  the  energy  eigenvalues  obtained  by 
diagonalizing  Eqn.  1  when  By  =  0  for  a  system  with  6 
spins.  The  ground  state  passes  through  three  level  cross¬ 
ings  as  the  strength  of  the  longitudinal  field  Bx  is  in¬ 
creased  from  Bx  =  0,  indicating  three  first-order  phase 
transitions  separating  four  distinct  spin  phases.  When 
Bx  is  set  at  a  phase  transition  and  By(t)  is  ramped  down, 
as  in  Fig.  1(b),  the  critical  gap  between  the  ground  and 
first  coupled  excited  state  closes  on  account  of  the  degen¬ 
eracy  of  spin  orderings  at  the  phase  transition. 

Long-range  interactions  allow  for  the  creation  of  many 
more  spin  phases  compared  with  a  nearest-neighbor-only 
Ising  model.  Consider  a  nearest-neighbor  AFM  model 
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FIG.  1.  (a)  Low-lying  energy  eigenvalues  of  Eqn.  1  for  By  =  0 
and  IV  =  6,  with  the  long-range  Jij  couplings  determined 
from  experimental  conditions  (see  text).  Level  crossings  (in¬ 
set)  indicate  the  presence  of  first-order  phase  transitions  in 
the  ground  state,  (b)  The  critical  gap  Ac  between  the  ground 
state  energy  Eg  and  first  coupled  excited  state  energy  shrinks 
to  zero  at  the  three  phase  transitions  (vertical  dashed  lines). 
Inset:  low-lying  energy  levels  of  Eqn.  1  as  the  transverse  field 
By  is  varied  with  Bx  =  0. 


with  N/2  “up”  spins,  N  total  spins,  and  a  ground-state 
ordering  |..  ■J-'N-tlt  ■■)x-  An  excited  state  at  Bx  =  0  may 
have  an  additional  spin  polarized  along  |4,)x,  either  by 
making  a  kink  of  type  |..  4-t44-t4-  ■■)x  or  a  spin  defect  of 
type  |..  4.f444-t  ..)x.  The  interaction  energy  cost  of  mak¬ 
ing  n  kinks  is  2 nJ ,  while  the  field  energy  gain  is  2 nBx. 
At  Bx/ J  =  1  multiple  energy  levels  intersect  to  give  a 
first-order  phase  transition.  Similarly,  the  energy  cost  of 
making  n  spin  defects  is  4 nJ  and  the  gain  is  2 nBXl  so  a 
second  phase  transition  occurs  at  Bx/J  =  2.  Only  three 
different  spin  phases  are  observable  as  Bx  is  varied  from 
0  — >  oo,  independent  of  IV,  and  there  is  a  large  degener¬ 
acy  of  spin  eigenstates  at  the  phase  transitions.  The  pres¬ 
ence  of  long-range  interactions  lifts  this  degeneracy  and 
admits  N/ 2  +  1  distinct  spin  phases,  with  {0, 1, ... ,  N/2} 
spins  in  state  |f)x- 

Observation  of  these  IV/ 2  +  1  classical  phases  is  made 
possible  only  by  quantum  fluctuations  in  our  quantum 
simulator.  The  second  term  in  Eq.  1  commutes  with 
the  first,  so  the  system  is  classical  when  By  =  0.  How¬ 
ever,  varying  Bx  across  a  phase  boundary  in  a  zero  tem¬ 
perature  system  does  not  lead  to  a  new  ground  state, 
since  thermal  fluctuations  are  required  to  drive  the  tran¬ 
sition.  Observation  of  spin  phases  at  intermediate  values 


of  Bx  and  zero  temperature  thus  requires  quantum  fluc¬ 
tuations  [30],  which  are  provided  by  the  non-commuting 
third  term  in  the  Hamiltonian  (Eqn.  1). 

The  effective  spin  system  is  encoded  in  a  linear  chain 
of  trapped  1(1Yb  ions  [31],  with  zero  effective  spin  tem¬ 
perature.  The  spin  states  1 1)^  and  |4.)z  are  represented 
by  the  hyperfine  clock  states  2S'1/2  | F  =  1  ,mp  =  0)  and 
|F  =  0,  rriF  =  0),  respectively,  which  have  a  frequency 
splitting  of  u>s/ 27T  =  12.642819  GHz  [32].  A  weak  mag¬ 
netic  field  of  ~  5  G  defines  the  quantization  axis.  The 
states  are  detected  by  illuminating  the  ions  with  laser 
light  resonant  with  the  2S'1/2  to  2Pi/2  cycling  transi¬ 
tion  at  369.5  nm  and  imaging  the  spin-dependent  flu¬ 
orescence.  Either  N  =  6  or  N  =  10  ions  are  confined 
in  a  three-layer  rf  Paul  trap  with  a  center-of-mass  axial 
trap  frequency  fz  =  0.7  MHz  and  transverse  frequencies 
fx  =  4.8  MHz  and  fy  =  4.6  MHz  and  interact  with  each 
other  via  their  collective  modes  of  motion. 

The  Ising  couplings  JZj  are  generated  by  globally  ir¬ 
radiating  the  ions  with  two  off-resonant  A  =  355  nm 
laser  beams  which  drive  stimulated  Raman  transitions 
[33,  34].  The  beams  intersect  at  right  angles  so  that  their 
wavevector  difference  Afc  points  along  the  cc-direction  of 
transverse  ion  motion,  perpendicular  to  the  linear  chain. 
Acousto-optic  modulators  imprint  beatnote  frequencies 
of  wg±/i  between  the  beams,  imparting  a  spin-dependent 
optical  dipole  force  at  frequency  p  [35] .  In  the  limit  where 
the  beatnotes  are  sufficiently  far  from  the  transverse  nor¬ 
mal  modes  com,  we  obtain  a  spin-spin  coupling  given  by 


Ij  '1  -  £1 ,  £1  n 


fc(Afc)2 
2  M 


Ebi,mbj,m 
M2  — 


(2) 


in  the  Lamb-Dicke  limit,  where  £2,;  is  the  Rabi  frequency 
of  the  ith  ion,  M  is  the  single  ion  mass,  and  6j>m  is  the 
normal-mode  transformation  matrix  for  ion  i  in  mode  m 
[17].  The  Ising  interactions  are  long-range  and  fall  off 
approximately  as  Jij  ~  \/\%  —  j\a,  with  a  =  0.94  for 
N  =  6  and  a  =  0.83  for  N  =  10. 

The  effective  transverse  and  longitudinal  magnetic 
fields  By(t)  and  Bx  drive  Rabi  oscillations  between  the 
spin  states  \\) z  and  |f)z.  Each  effective  field  is  gener¬ 
ated  by  a  pair  of  Raman  laser  beams  with  a  beatnote 
frequency  of  cos,  with  the  field  amplitude  determined  by 
the  beam  intensities.  The  field  directions  are  controlled 
through  the  beam  phases  relative  to  the  average  phase  ip 
of  the  two  sidebands  which  give  rise  to  the  axax  inter¬ 
action  in  Eqn.  2.  In  particular,  an  effective  field  phase 
offset  of  0°  (90°)  relative  to  <p  generates  a  ay  (ax)  inter¬ 
action. 

Each  experiment  begins  with  3  ms  of  Doppler  cool¬ 
ing,  followed  by  optical  pumping  to  the  state  |444-  ..) z 
and  100  /zs  of  Raman  sideband  cooling  that  prepares  the 
motion  of  all  modes  along  Afc  in  the  Lamb-Dicke  limit. 
The  spins  are  then  coherently  rotated  into  the  equatorial 
plane  of  the  Bloch  sphere  so  that  they  point  along  the 
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FIG.  2.  (Color)  (a)  Magnetization  (mz  =  Nf  —  N±)  of  a  chain 
of  6  ions  for  increasing  axial  field  strength.  Red,  solid:  pre¬ 
dicted  magnetization  value,  with  the  step  locations  indicat¬ 
ing  the  first-order  phase  transitions.  Blue  diamonds:  average 
magnetization  of  4000  experiments  performed  at  each  value 
of  Bx.  Black,  dashed:  magnetization  of  the  most  probable 
state  (see  inset)  found  at  each  Bx  value,  (b)  Camera  images 
of  the  ground  states  found  at  each  step  in  (a):  [414141") ,  and 

limn),  (m,  =  o),  imm>,  k  =  -2),  mtm>,  and 

|illtii)z  (m,  =  -4),  and  |444444)_  K  =  -6). 

direction  B  =  Bxx  +  Bv(0)y,  where  Bx  is  varied  between 
different  simulations.  The  Hamiltonian  (Eqn.  1)  is  then 
switched  on  at  t  =  0  with  the  chosen  value  of  Bx  and 
By( 0)  =  5Jmax,  where  Jmax  is  the  largest  spin-spin  cou¬ 
pling  (typically  27tx0.6-0.7  kHz).  The  transverse  field  is 
ramped  down  to  By  ss  0  exponentially  with  a  time  con¬ 
stant  of  600  fj, s  and  a  total  time  of  3  ms,  with  the  ramp 
times  chosen  to  minimize  decoherence  effects  while  max¬ 
imizing  adiabaticity.  At  t  =  3  ms,  the  Hamiltonian  is 
switched  off  and  the  a;— component  of  each  spin  is  mea¬ 
sured  by  applying  a  global  tt/2  rotation  about  the  y  axis, 
illuminating  the  ions  with  resonant  light,  and  imaging 
the  spin-dependent  fluorescence  using  an  intensified  CCD 
camera  (Princeton  Instruments  PIMax3:1024i).  The  ex¬ 
periments  are  repeated  4000  times  to  determine  the  prob¬ 
ability  of  each  possible  spin  configuration. 

The  camera  measurement  cycle  for  each  experiment 


takes  10  ms,  with  3  ms  of  fluorescence  collection,  3  ms  of 
phosphor  decay  from  the  CCD  intensifier,  and  3.6  ms  of 
readout  time.  For  all  data,  we  restrict  the  region  of  inter¬ 
est  on  the  CCD  to  an  approximately  128  x  7  bin  window 
to  speed  readout  from  the  device.  We  compensate  for 
detection  errors  (e  =  7%  for  a  single  spin)  by  multiplying 
a  matrix  describing  the  expected  multi-spin  error  by  the 
vector  containing  the  measured  probability  of  each  spin 
configuration  [36]. 

We  investigate  the  order  parameter  of  net  magnetiza¬ 
tion,  mz  =  Nf  —  JVj.  as  we  vary  the  strength  of  the  longi¬ 
tudinal  field  Bx.  The  magnetization  of  the  ground  state 
spin  ordering  of  Eq.  1  is  expected  to  yield  a  staircase 
with  sharp  steps  at  the  phase  transitions  (red  line  in  Fig. 
2(a))  when  By  =  0  [29].  The  experimental  data  (blue 
points  in  Fig.  2(a))  show  an  averaged  magnetization  with 
heavily  broadened  steps  due  largely  to  non-adiabatic  evo¬ 
lution  during  the  exponential  ramp  that  populate  excited 
states.  However,  the  ground  state  spin  configuration  at 
each  value  of  Bx  may  be  extracted  by  looking  at  the 
probability  distribution  of  all  spin  states  and  selecting 
the  most  prevalent  state  (inset  of  Fig.  2(a))  [37].  The 
magnetization  of  the  spin  states  found  by  this  method 
(black  points  in  Fig.  2(a))  agree  well  with  the  theoreti¬ 
cal  prediction. 

Direct  images  of  the  ground  state  spin  phases  are 
shown  in  Fig.  2(b),  with  approximately  1000  averaged 
images  per  line.  Each  box  contains  an  ion  that  scatters 
many  photons  when  in  the  state  |f),  and  essentially  no 
photons  when  in  the  state  |4)z.  For  magnetizations  of  0 
and  —4,  two  ground  state  orderings  are  observed  due  to 
the  left-right  symmetry  of  the  spin-spin  interactions.  In 
cases  where  left-right  reflection  of  the  spin  chain  produces 
a  new  spin  ordering  with  degenerate  energy,  the  summed 
probability  of  both  states  is  used  when  determining  the 
most  probable  state. 


B*/  ax 


FIG.  3.  (Color  Online)  Probabilities  of  the  4  different  ground 
state  spin  phases  as  Bx  is  varied  in  a  6-ion  system.  Black 

dots:  1414141")  z  or  1 141414) z-  Red  squares:  [414414),  ■  Blue 
diamonds:  1 444444) z  or  1 444444) z-  Green  triangles:  1444444), 

The  probability  of  populating  the  four  distinct  ground 
states  as  Bx  is  varied,  shown  in  Fig.  3,  further  probes 
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FIG.  4.  (Color)  (a)  Magnetization  of  a  chain  of  10  ions  for  in¬ 
creasing  axial  field  strength.  The  red,  blue,  and  black  curves 
correspond  to  the  theoretical  magnetization,  average  mea¬ 
sured  magnetization,  and  magnetization  of  the  most  probable 
state  (respectively)  for  increasing  Bx.  (b)  Camera  images  of 
the  theoretical  ground  states  found  experimentally  at  each 
step  in  (a). 


the  multiple  phase  transitions.  Fig.  1  shows  that  as 
Bx/ Jm ax  approaches  a  phase  transition,  the  energy  lev¬ 
els  of  the  spin  eigenstates  on  opposing  sides  of  the  transi¬ 
tion  cross,  and  the  critical  gap  Ac  between  them  shrinks. 
The  smooth  crossover  between  the  four  distinct  ground- 
state  spin  phases  in  Fig.  3  therefore  arises  from  increas¬ 
ingly  diabatic  simulations  near  the  energy-level  crossings 
at  the  first-order  transitions.  Comparing  with  Fig.  1(b), 
we  find  higher  probabilities  of  creating  the  theoretical 
ground  state  when  the  critical  gap  Ac  is  large.  Simula¬ 
tions  that  are  more  adiabatic  could  be  expected  to  give 
both  larger  maximum  probabilities  in  each  of  the  four 
different  ground  states  and  sharper  crossovers  between 
the  state  populations  near  the  phase  transitions. 

Figure  4  plots  the  magnetization  of  a  10- ion  chain  for 
increasing  strengths  of  the  longitudinal  field  Bx,  along 
with  the  associated  camera  images  of  the  N/2  +  1  =  6 
distinct  ground  state  phases.  Similarly  to  the  6-ion  case, 
we  experimentally  determine  the  ground  state  ordering 
by  selecting  the  most  prevalent  of  the  210  =  1024  possible 


states.  Because  the  critical  gap  Ac  is  much  smaller  for  a 
10-ion  system,  the  reduction  of  the  transverse  field  By(t) 
is  much  less  adiabatic,  and  only  ~  5%  of  the  total  pop¬ 
ulation  lands  in  the  theoretical  ground  state  for  much  of 
the  data  in  Fig.  4(a).  We  are  therefore  more  sensitive  to 
experimental  error  sources  (of  order  a  few  percent)  such 
as  imperfect  initialization  and  slow  drifts  in  the  strength 
of  Ji  j,  as  well  as  statistically-limited  quantum  projection 
noise  of  order  «  1%. 

For  the  10-ion  chain,  Fig.  4(b)  shows  the  interesting 
spin  structure  that  emerges  as  Bx  is  varied.  For  a  given 
Bx  and  associated  number  of  bright  ions  q ,  the  ground- 
state  spin  configuration  of  Eqn.  1  (with  By  =  0)  solves 
the  problem  of  finding  the  lowest  energy  arrangement 
of  q  charged  particles  on  N  lattice  sites.  The  creation 
of  such  periodic  spin  structures  realizes  a  generalized 
Wigner  crystal,  which  is  the  configuration  adopted  by 
a  low-density  electron  gas  when  the  Coulomb  repulsion 
dominates  over  the  kinetic  energy  [38].  As  the  system 
size  N  — >  oo  and  Bx  is  increased,  the  magnetization  re¬ 
veals  a  fractal  staircase  structure  that  arises  since  every 
rational  filling  factor  (of  which  there  are  infinitely  many) 
is  the  ground  state  for  some  value  of  Bx  [29]. 

In  conclusion,  we  have  performed  a  quantum  simula¬ 
tion  of  an  Ising  model  with  both  transverse  and  longi¬ 
tudinal  magnetic  fields  within  a  linear  chain  of  trapped 
ions.  Long-range  spin-spin  couplings  enabled  the  exis¬ 
tence  of  a  large  number  of  ground  state  spin  phases  at 
zero  temperature,  while  quantum  fluctuations  inherent  in 
our  quantum  simulations  enabled  their  observation.  As 
the  system  size  is  scaled  up,  such  quantum  simulations 
can  begin  to  solve  diverse  problems  in  many-body  physics 
where  classical  calculation  becomes  intractable. 
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